Let us set some global options for all code chunks in this document.

# Set seed for reproducibility
set.seed(1982) 
# Set global options for all code chunks
knitr::opts_chunk$set(
  # Disable messages printed by R code chunks
  message = TRUE,    
  # Disable warnings printed by R code chunks
  warning = TRUE,    
  # Show R code within code chunks in output
  echo = TRUE,        
  # Include both R code and its results in output
  include = TRUE,     
  # Evaluate R code chunks
  eval = TRUE,       
  # Enable caching of R code chunks for faster rendering
  cache = FALSE,      
  # Align figures in the center of the output
  fig.align = "center",
  # Enable retina display for high-resolution figures
  retina = 2,
  # Show errors in the output instead of stopping rendering
  error = TRUE,
  # Do not collapse code and output into a single block
  collapse = FALSE
)
# Start the figure counter
fig_count <- 0
# Define the captioner function
captioner <- function(caption) {
  fig_count <<- fig_count + 1
  paste0("Figure ", fig_count, ": ", caption)
}
# Define the function to truncate a number to two decimal places
truncate_to_two <- function(x) {
  floor(x * 100) / 100
}
my.get.roots <- function(order, beta) {
  mt <- get(paste0("m", order, "table"))
  rb <- rep(0, order + 1)
  rc <- rep(0, order)
  if(order == 1) {
      rc = approx(mt$beta, mt[[paste0("rc")]], beta)$y
    } else {
      rc = sapply(1:order, function(i) {
        approx(mt$beta, mt[[paste0("rc.", i)]], beta)$y
      })
    }
    rb = sapply(1:(order+1), function(i) {
      approx(mt$beta, mt[[paste0("rb.", i)]], xout = beta)$y
    })
    factor = approx(mt$beta, mt$factor, xout = beta)$y
  return(list(rb = rb, rc = rc, factor = factor))
}

m1table <- rSPDE:::m1table
m2table <- rSPDE:::m2table
m3table <- rSPDE:::m3table
m4table <- rSPDE:::m4table

poly_from_roots <- function(roots) {
  coef <- 1
  for (r in roots) {
    coef <- convolve(coef, c(1, -r), type = "open")
  }
  return(coef) # returns in increasing order like a+bx+cx^2+...
}
 
# Group complex roots into conjugate pairs
get_conjugate_pairs <- function(roots) {
  used <- rep(FALSE, length(roots))
  pairs <- list()
  for (i in seq_along(roots)) {
    if (!used[i]) {
      conj_root <- Conj(roots[i])
      # Find its conjugate (within tolerance)
      j <- which(!used & abs(Re(roots) - Re(conj_root)) < 1e-8 & abs(Im(roots) - Im(conj_root)) < 1e-8)
      j <- j[j != i]
      if (length(j) > 0) {
        used[c(i, j[1])] <- TRUE
        pairs[[length(pairs) + 1]] <- c(roots[i], roots[j[1]])
      }
    }
  }
  return(pairs)
}

compute_sum_poly <- function(factor, pr_roots, pl_roots, cte) {
  
  pr_coef <- c(0, poly_from_roots(pr_roots)) # in decreasing order like x^n+bx^(n-1)+cx^(n-2)+...
  pl_coef <- poly_from_roots(pl_roots) # in decreasing order like x^n+bx^(n-1)+cx^(n-2)+...
  pr_plus_pl_coef <- factor * pr_coef + cte * pl_coef
  return(pr_plus_pl_coef)
}

compute_real_roots_and_complex_coef <- function(pr_plus_pl_coef) {

  pr_plus_pl_roots <- polyroot(rev(pr_plus_pl_coef))
  
  real_roots <- Re(pr_plus_pl_roots[abs(Im(pr_plus_pl_roots)) < 1e-8])
  complex_roots <- pr_plus_pl_roots[abs(Im(pr_plus_pl_roots)) >= 1e-8]
  
  complex_poly_coefs <- list()
  if (length(complex_roots) > 0) {
    pairs <- get_conjugate_pairs(complex_roots)
    for (pair in pairs) {
      coef <- rev(Re(poly_from_roots(c(pair[1], pair[2]))))  # returns in x^2 + bx + c order
      complex_poly_coefs[[length(complex_poly_coefs) + 1]] <- round(coef, 12)
    }
  } else {
    complex_poly_coefs <- list()  # No complex roots
  }
  
  return(list(new_factor = pr_plus_pl_coef[1], 
              pr_plus_pl_roots = pr_plus_pl_roots, 
              real_roots = real_roots, 
              complex_poly_coefs = complex_poly_coefs))
}

my.fractional.operators <- function(L, # kappa^2C + G
                                 beta,
                                 C,
                                 scale.factor, # kappa^2
                                 m = 1,
                                 time_step) {

  C <- Matrix::Diagonal(dim(C)[1], rowSums(C)) # lumped
  Ci <- Matrix::Diagonal(dim(C)[1], 1 / rowSums(C)) # lumped 
  I <- Matrix::Diagonal(dim(C)[1])
  L <- L / scale.factor # C + G/kappa^2
  LCi <- L %*% Ci
  roots <- my.get.roots(m, beta)
  Pl.roots <- roots$rb
  Pr.roots <- roots$rc
  factor <- roots$factor
  new_factor_and_roots <- compute_real_roots_and_complex_coef(compute_sum_poly(factor, Pr.roots, Pl.roots, time_step))
  
  new_real_roots <- new_factor_and_roots$real_roots
  new_factor <- new_factor_and_roots$new_factor
  complex_poly_coefs <- new_factor_and_roots$complex_poly_coefs
  
  Pr_plus_Pl.factors <- list()
  if (length(new_real_roots) >= 1) {
    for (i in 1:length(new_real_roots)) {
      Pr_plus_Pl.factors[[i]] <- LCi - new_real_roots[i] * I
    }
  }
  length_now <- length(Pr_plus_Pl.factors)
  if(length(complex_poly_coefs) >= 1) {
    for (i in 1:length(complex_poly_coefs)) {
    }
    Pr_plus_Pl.factors[[length_now + i]] <- LCi %*% LCi + complex_poly_coefs[[i]][2] * LCi + complex_poly_coefs[[i]][3] * I
  }
  
  Pr.factors <- list()
  Pr.factors[[1]] <- I - LCi * roots$rc[1]

  if (length(roots$rc) > 1) {
    for (i in 2:length(roots$rc)) {
      Pr.factors[[i]] <- I - LCi * roots$rc[i]
    }
  }
  
  #Pr.factors[[length(Pr.factors) + 1]] <- I
   
  output <- list(
    Ci = Ci,
    C = C,
    LCi = LCi,
    L = L,
    m = m,
    beta = beta,
    factor = factor,
    new_factor = new_factor,
    Pr.factors = Pr.factors,
    Pr_plus_Pl.factors = Pr_plus_Pl.factors
  )
  return(output)
}
# install.packages("INLA",repos=c(getOption("repos"),INLA="https://inla.r-inla-download.org/R/testing"), dep=TRUE)
# inla.upgrade(testing = TRUE)
# remotes::install_github("inlabru-org/inlabru", ref = "devel")
# remotes::install_github("davidbolin/rspde", ref = "devel")
# remotes::install_github("davidbolin/metricgraph", ref = "devel")
library(INLA)
## Loading required package: Matrix
## This is INLA_25.05.01 built 2025-05-01 18:43:33 UTC.
##  - See www.r-inla.org/contact-us for how to get help.
##  - List available models/likelihoods/etc with inla.list.models()
##  - Use inla.doc(<NAME>) to access documentation
##  - Consider upgrading R-INLA to testing[25.05.04-1] or stable[24.12.11].
library(inlabru)
## Loading required package: fmesher
library(rSPDE)
## This is rSPDE 2.5.1
## - See https://davidbolin.github.io/rSPDE for vignettes and manuals.
library(MetricGraph)
## This is MetricGraph 1.4.1
## - See https://davidbolin.github.io/MetricGraph for vignettes and manuals.
## 
## Attaching package: 'MetricGraph'
## The following object is masked from 'package:stats':
## 
##     filter
library(grateful)

library(plotly)
## Loading required package: ggplot2
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout

We want to solve the fractional diffusion equation \[\begin{equation} \label{eq:maineq} \partial_t u+(\kappa^2-\Delta_\Gamma)^{\frac{\alpha}{2}} u=f \text { on } \Gamma \times(0, T), \quad u(0)=u_0 \text { on } \Gamma, \end{equation}\] where \(u\) satisfies the Kirchhoff vertex conditions \[\begin{equation} \label{eq:Kcond} \left\{\phi\in C(\Gamma)\;\Big|\; \forall v\in V: \sum_{e\in\mathcal{E}_v}\partial_e \phi(v)=0 \right\} \end{equation}\]

If \(f=0\), then the solution is given by \[\begin{equation} \label{eq:sol_reprentation} u(s,t) = \displaystyle\sum_{j\in\mathbb{N}}e^{-\lambda^{\frac{\alpha}{2}}_jt}\left(u_0, e_j\right)_{L_2(\Gamma)}e_j(s). \end{equation}\]

# Function to build a tadpole graph and create a mesh
gets_graph_tadpole <- function(h){
  edge1 <- rbind(c(0,0),c(1,0))
  theta <- seq(from=-pi,to=pi,length.out = 100)
  edge2 <- cbind(1+1/pi+cos(theta)/pi,sin(theta)/pi)
  edges = list(edge1, edge2)
  graph <- metric_graph$new(edges = edges)
  graph$build_mesh(h = h)
  return(graph)
}

Let \(\Gamma_T = (\mathcal{V},\mathcal{E})\) characterize the tadpole graph with \(\mathcal{V}= \{v_1,v_2\}\) and \(\mathcal{E}= \{e_1,e_2\}\). The left edge \(e_1\) has length 1 and the circular edge \(e_2\) has length 2. As discussed before, a point on \(e_1\) is parameterized via \(s=\left(e_1, t\right)\) for \(t \in[0,1]\) and a point on \(e_2\) via \(s=\left(e_2, t\right)\) for \(t\in[0,2]\). One can verify that \(-\Delta_\Gamma\) has eigenvalues \(0,\left\{(i \pi / 2)^2\right\}_{i \in \mathbb{N}}\) and \(\left\{(i \pi / 2)^2\right\}_{2 i \in \mathbb{N}}\) with corresponding eigenfunctions \(\phi_0\), \(\left\{\phi_i\right\}_{i \in \mathbb{N}}\), and \(\left\{\psi_i\right\}_{2 i \in \mathbb{N}}\) given by \(\phi_0(s)=1 / \sqrt{3}\) and \[\begin{equation*} \phi_i(s)=C_{\phi, i}\begin{cases} -2 \sin (\frac{i\pi}{2}) \cos (\frac{i \pi t}{2}), & s \in e_1, \\ \sin (i \pi t / 2), & s \in e_2, \end{cases}, \quad \psi_i(s)=\frac{\sqrt{3}}{\sqrt{2}} \begin{cases} (-1)^{i / 2} \cos (\frac{i \pi t}{2}), & s \in e_1, \\ \cos (\frac{i \pi t}{2}), & s \in e_2, \end{cases}, \end{equation*}\] where \(C_{\phi, i}=1\) if \(i\) is even and \(C_{\phi, i}=1 / \sqrt{3}\) otherwise. Moreover, these functions form an orthonormal basis for \(L_2(\Gamma_T)\).

# Function to compute the eigenfunctions 
tadpole.eig <- function(k,graph){
x1 <- c(0,graph$get_edge_lengths()[1]*graph$mesh$PtE[graph$mesh$PtE[,1]==1,2]) 
x2 <- c(0,graph$get_edge_lengths()[2]*graph$mesh$PtE[graph$mesh$PtE[,1]==2,2]) 

if(k==0){ 
  f.e1 <- rep(1,length(x1)) 
  f.e2 <- rep(1,length(x2)) 
  f1 = c(f.e1[1],f.e2[1],f.e1[-1], f.e2[-1]) 
  f = list(phi=f1/sqrt(3)) 
  
} else {
  f.e1 <- -2*sin(pi*k*1/2)*cos(pi*k*x1/2) 
  f.e2 <- sin(pi*k*x2/2)                  
  
  f1 = c(f.e1[1],f.e2[1],f.e1[-1], f.e2[-1]) 
  
  if((k %% 2)==1){ 
    f = list(phi=f1/sqrt(3)) 
  } else { 
    f.e1 <- (-1)^{k/2}*cos(pi*k*x1/2)
    f.e2 <- cos(pi*k*x2/2)
    f2 = c(f.e1[1],f.e2[1],f.e1[-1],f.e2[-1]) 
    f <- list(phi=f1,psi=f2/sqrt(3/2))
  }
}

return(f)
}

Implementation of \(u\)

h <- 0.01
graph <- gets_graph_tadpole(h = h)
## Starting graph creation...
## LongLat is set to FALSE
## Creating edges...
## Setting edge weights...
## Computing bounding box...
## Setting up edges
## Merging close vertices
## Total construction time: 0.16 secs
## Creating and updating vertices...
## Storing the initial graph...
## Computing the relative positions of the edges...
T_final <- 2
time_step <- 0.01
time_seq <- seq(0, T_final, by = time_step)
# Compute the FEM matrices
graph$compute_fem()
G <- graph$mesh$G
C <- graph$mesh$C
I <- Matrix::Diagonal(nrow(C))
x <- graph$mesh$V[, 1]
y <- graph$mesh$V[, 2]
edge_number <- graph$mesh$VtE[, 1]
pos <- sum(edge_number == 1)+1
order_to_plot <- function(v)return(c(v[1], v[3:pos], v[2], v[(pos+1):length(v)], v[2]))
weights <- graph$mesh$weights
# Initial condition
# U_0 <- 10*exp(-((x-1)^2 + (y)^2))
kappa <- 1
alpha <- 0.5 # from 0.5 to 2
m = 4
beta <- alpha/2
L <- kappa^2*C + G
my_op <- my.fractional.operators(L, beta, C, scale.factor = kappa^2, m = m, time_step)

my.solver <- function(obj, v){
  Pr.factors <- obj$Pr.factors
  Pr_plus_Pl.factors <- obj$Pr_plus_Pl.factors
  m <- obj$m
  C <- obj$C
  Ci <- obj$Ci
  factor <- obj$factor
  new_factor <- obj$new_factor
  if (m==1){
    temp <- solve(Pr_plus_Pl.factors[[2]], 
                  Pr.factors[[1]] %*% solve(
                    Pr_plus_Pl.factors[[1]], 
                    C%*%v))
  } else if (m==2){
    temp <- solve(Pr_plus_Pl.factors[[3]], 
                  Pr.factors[[2]] %*% solve(
                    Pr_plus_Pl.factors[[2]], 
                    Pr.factors[[1]] %*% solve(
                      Pr_plus_Pl.factors[[1]], 
                      C%*%v)))
  } else if (m==3){
    temp <- solve(Pr_plus_Pl.factors[[4]], 
                  Pr.factors[[3]] %*% solve(
                    Pr_plus_Pl.factors[[3]], 
                    Pr.factors[[2]] %*% solve(
                      Pr_plus_Pl.factors[[2]], 
                      Pr.factors[[1]] %*% solve(
                        Pr_plus_Pl.factors[[1]], 
                        C%*%v))))
  } else if (m==4){
    temp <- solve(Pr_plus_Pl.factors[[5]], 
                  Pr.factors[[4]] %*% solve(
                    Pr_plus_Pl.factors[[4]], 
                    Pr.factors[[3]] %*% solve(
                      Pr_plus_Pl.factors[[3]], 
                      Pr.factors[[2]] %*% solve(
                        Pr_plus_Pl.factors[[2]], 
                        Pr.factors[[1]] %*% solve(
                          Pr_plus_Pl.factors[[1]], 
                          C%*%v)))))
  }
  return((factor/new_factor) * Ci %*% temp)
}
op <- fractional.operators(L, beta, C, scale.factor = kappa^2, m = m)
Pl <- op$Pl
Pr <- op$Pr
Ci <- op$Ci
C <- op$C
# Parameters to construct U_0
N_finite <- 4 # choose an even number
adjusted_N_finite <- N_finite + N_finite/2 + 1
EIGENVAL_ALPHA <- NULL
EIGENFUN <- NULL       
for (j in 0:N_finite) {
    lambda_j_alpha <- (kappa^2 + (j*pi/2)^2)^(alpha/2)
    e_j <- tadpole.eig(j,graph)$phi
    EIGENVAL_ALPHA <- c(EIGENVAL_ALPHA, lambda_j_alpha)         
    EIGENFUN <- cbind(EIGENFUN, e_j)            
    if (j>0 && (j %% 2 == 0)) {
      lambda_j_alpha <- (kappa^2 + (j*pi/2)^2)^(alpha/2)
      e_j <- tadpole.eig(j,graph)$psi
      EIGENVAL_ALPHA <- c(EIGENVAL_ALPHA, lambda_j_alpha)         
      EIGENFUN <- cbind(EIGENFUN, e_j)           
    }
}

# Building the initial condition as \sum coeff_j EIGENFUN_j
coeff <- 1*(1:adjusted_N_finite)^-1
coeff[-5] <- 0
# lower_zeroer <- 10  # choose an even number
# adjusted_lower_zeroer <- lower_zeroer + lower_zeroer/2 + 1
# coeff[adjusted_lower_zeroer:adjusted_N_finite] <- 0
U_0 <- EIGENFUN %*% coeff

# Building the true solution as \sum coeff_j EIGENFUN_j e^{-\lambda_j^{\frac{\alpha}{2}}t}
U_true <- matrix(NA, nrow = length(x), ncol = length(time_seq))
for (k in 1:length(time_seq)) {
  aux_k <- rep(0, length(x))
  for (j in 1:adjusted_N_finite) {
    aux_k <- aux_k + exp(-time_seq[k]*EIGENVAL_ALPHA[j])*coeff[j]*EIGENFUN[, j]
  }
  U_true[, k] <- aux_k
}
U_approx <- matrix(NA, nrow = nrow(C), ncol = length(time_seq))
U_approx[, 1] <- U_0

# Time-stepping loop
for (k in 1:(length(time_seq) - 1)) {
  U_approx[, k + 1] <- as.matrix(my.solver(my_op, U_approx[, k]))
}
{r}
Pr.apply.mult <- function(v){return(Pr.mult(op, v))}
Pl.apply.solve <- function(v){return(Pl.solve(op, v))}
PliC <- apply(C, 2, Pl.apply.solve) # PlC^-1
# Precompute the LHS1 matrix
aux <- apply(PliC, 2, Pr.apply.mult) #PrPl^-1C
LHS <- aux + time_step * Matrix::Diagonal(nrow(C)) 

# Initialize U matrix to store solution at each time step
U_approx <- matrix(NA, nrow = nrow(C), ncol = length(time_seq))
U_approx[, 1] <- U_0

# Time-stepping loop
for (k in 1:(length(time_seq) - 1)) {
  # Compute the right-hand side for the second equation
  RHS <- aux %*% U_approx[, k]
  U_approx[, k + 1] <- as.matrix(solve(LHS, RHS))
}

We arrive at the scheme

\[\begin{equation} (P_r^TC+\tau P_\ell^T)U^{k+1} = P_r^TCU^k \label{eq:scheme} \end{equation}\] where \[\begin{equation} P_r = c_m\prod_{i=1}^m (I-r_{1i}C^{-1}L)\quad\text{and}\quad P_\ell = b_{m+1}C\prod_{j=1}^{m+1} (I-r_{2j}C^{-1}L) \end{equation}\] Observe that \[\begin{equation} P_r^T = c_m\prod_{i=m}^1 (I-r_{1i}LC^{-1})\quad\text{and}\quad P_\ell^T = b_{m+1}\prod_{j=m+1}^{1} (I-r_{2j}LC^{-1})\cdot C \end{equation}\] Replacing these two in our scheme we get \[\begin{equation} \left(c_m\prod_{i=m}^1 (I-r_{1i}LC^{-1})+\tau b_{m+1}\prod_{j=m+1}^{1} (I-r_{2j}LC^{-1})\right)CU^{k+1} = c_m\prod_{i=m}^1 (I-r_{1i}LC^{-1})\cdot CU^k \end{equation}\] We can equivalently write this as \[\begin{equation} \left(\dfrac{c_m}{b_{m+1}}\prod_{i=1}^m (I-r_{1i}LC^{-1})+\tau \prod_{j=1}^{m+1} (I-r_{2j}LC^{-1})\right)CU^{k+1} = \dfrac{c_m}{b_{m+1}}\prod_{i=1}^m (I-r_{1i}LC^{-1})\cdot CU^k \end{equation}\]

We need to compute the roots of the polynomial \[\begin{equation} R(x) = \dfrac{c_m}{b_{m+1}} \prod_{i=1}^m (x-r_{1i})+\tau \prod_{j=1}^{m+1} (x-r_{2j}) \end{equation}\] Say \(R\) has leading coefficient \(\tau\) and \(m+1\) roots \(R_k\) for \(k=1,\ldots,m+1\). That is, \[\begin{equation} R(x) = \tau\prod_{k=1}^{m+1} (x-R_k) \end{equation}\] We can then write our scheme as \[\begin{equation} \left(\tau\prod_{k=1}^{m+1} (I-R_kLC^{-1})\right)CU^{k+1} = \dfrac{c_m}{b_{m+1}}\prod_{i=1}^m (I-r_{1i}LC^{-1})\cdot CU^k \end{equation}\] That is, \[\begin{equation} U^{k+1} = \dfrac{c_m}{b_{m+1}}\dfrac{1}{\tau}C^{-1}\prod_{k=1}^{m+1} (I-R_kLC^{-1})^{-1}\prod_{i=1}^m (I-r_{1i}LC^{-1})\cdot CU^k \end{equation}\]

What if \[\begin{equation} R(x) = \tau(x^2+ax+b)\prod_{k=1}^{m-1} (x-R_k) \end{equation}\] then

{r}
LHS <- t(Pr)%*% C + time_step * t(Pl)
# Initialize U matrix to store solution at each time step
U_approx <- matrix(NA, nrow = nrow(C), ncol = length(time_seq))
U_approx[, 1] <- U_0

# Time-stepping loop
for (k in 1:(length(time_seq) - 1)) {
  # Compute the right-hand side for the second equation
  RHS <- t(Pr)%*% C %*% U_approx[, k]
  U_approx[, k + 1] <- as.matrix(solve(LHS, RHS))
}
x <- order_to_plot(x)
y <- order_to_plot(y)
max_error_at_each_time <- apply((U_true - U_approx)^2, 2, mean)

U_true <- apply(U_true, 2, order_to_plot)
U_approx <- apply(U_approx, 2, order_to_plot)



# Create interactive plot
fig <- plot_ly()


# Add second line (max_error_at_each_time)
fig <- fig %>% add_trace(
  x = ~time_seq, y = ~max_error_at_each_time, type = 'scatter', mode = 'lines+markers',
  line = list(color = 'red', width = 2),
  marker = list(size = 4),
  name = "Max Error True and Approx 2"
)

# Layout
fig <- fig %>% layout(
  title = "Max Error at Each Time Step",
  xaxis = list(title = "Time"),
  yaxis = list(title = "Max Error"),
  legend = list(x = 0.1, y = 0.9)
)


plot_data <- data.frame(
  x = rep(x, times = ncol(U_true)),
  y = rep(y, times = ncol(U_true)),
  z_true = as.vector(U_true),
  z_approx = as.vector(U_approx),
  frame = rep(time_seq, each = length(x))
)

# Compute axis limits
x_range <- range(x)
y_range <- range(y)
z_range <- range(c(U_true, U_approx))

# Initial plot setup (first frame only)
p <- plot_ly(plot_data, frame = ~frame) %>%
  add_trace(
    x = ~x, y = ~y, z = ~z_true,
    type = "scatter3d", mode = "lines",
    name = "True",
    line = list(color = "blue", width = 2)
  ) %>%
  add_trace(
    x = ~x, y = ~y, z = ~z_approx,
    type = "scatter3d", mode = "lines",
    name = "Approx",
    line = list(color = "red", width = 2)
  ) %>%
  layout(
    scene = list(
      xaxis = list(title = "x", range = x_range),
      yaxis = list(title = "y", range = y_range),
      zaxis = list(title = "Value", range = z_range),
      aspectratio = list(x = 2.4, y = 1.2, z = 1.2),
           camera = list(
      eye = list(x = 1.5, y = 1.5, z = 1),  # Adjust the viewpoint
      center = list(x = 0, y = 0, z = 0))
    ),
    updatemenus = list(
      list(
        type = "buttons", showactive = FALSE,
        buttons = list(
          list(label = "Play", method = "animate",
               args = list(NULL, list(frame = list(duration = 100, redraw = TRUE), fromcurrent = TRUE))),
          list(label = "Pause", method = "animate",
               args = list(NULL, list(mode = "immediate", frame = list(duration = 0), redraw = FALSE)))
        )
      )
    ),
    title = "Time: 0"
  )

# Convert to plotly object with frame info
pb <- plotly_build(p)

# Inject custom titles into each frame
for (i in seq_along(pb$x$frames)) {
  t <- time_seq[i]
  err <- signif(max_error_at_each_time[i], 4)
  pb$x$frames[[i]]$layout <- list(title = paste0("Time: ", t, " | Max Error: ", err))
}
fig  # Display the plot
pb

Figure 1: Caption

LS0tCnRpdGxlOiAiU29sdmluZyBhIHBhcmFib2xpYyBlcXVhdGlvbiIKZGF0ZTogIkNyZWF0ZWQ6IDIwLTA0LTIwMjUuIExhc3QgbW9kaWZpZWQ6IGByIGZvcm1hdChTeXMudGltZSgpLCAnJWQtJW0tJVkuJylgIgpvdXRwdXQ6CiAgaHRtbF9kb2N1bWVudDoKICAgIG1hdGhqYXg6ICJodHRwczovL2Nkbi5qc2RlbGl2ci5uZXQvbnBtL21hdGhqYXhAMy9lczUvdGV4LW1tbC1jaHRtbC5qcyIKICAgIGhpZ2hsaWdodDogcHlnbWVudHMKICAgIHRoZW1lOiBmbGF0bHkKICAgIGNvZGVfZm9sZGluZzogc2hvdyAjIGNsYXNzLnNvdXJjZSA9ICJmb2xkLWhpZGUiIHRvIGhpZGUgY29kZSBhbmQgYWRkIGEgYnV0dG9uIHRvIHNob3cgaXQKICAgIGRmX3ByaW50OiBwYWdlZAogICAgIyB0b2M6IHRydWUKICAgICMgdG9jX2Zsb2F0OgogICAgIyAgIGNvbGxhcHNlZDogdHJ1ZQogICAgIyAgIHNtb290aF9zY3JvbGw6IHRydWUKICAgIG51bWJlcl9zZWN0aW9uczogZmFsc2UKICAgIGZpZ19jYXB0aW9uOiB0cnVlCiAgICBjb2RlX2Rvd25sb2FkOiB0cnVlCmFsd2F5c19hbGxvd19odG1sOiB0cnVlCmJpYmxpb2dyYXBoeTogCiAgLSByZWZlcmVuY2VzLmJpYgogIC0gZ3JhdGVmdWwtcmVmcy5iaWIKaGVhZGVyLWluY2x1ZGVzOgogIC0gXG5ld2NvbW1hbmR7XGFyfXtcbWF0aGJie1J9fQogIC0gXG5ld2NvbW1hbmR7XGxsYXZ9WzFde1xsZWZ0XHsjMVxyaWdodFx9fQogIC0gXG5ld2NvbW1hbmR7XHBhcmV9WzFde1xsZWZ0KCMxXHJpZ2h0KX0KICAtIFxuZXdjb21tYW5ke1xOY2FsfXtcbWF0aGNhbHtOfX0KICAtIFxuZXdjb21tYW5ke1xWY2FsfXtcbWF0aGNhbHtWfX0KICAtIFxuZXdjb21tYW5ke1xFY2FsfXtcbWF0aGNhbHtFfX0KICAtIFxuZXdjb21tYW5ke1xXY2FsfXtcbWF0aGNhbHtXfX0KLS0tCgpgYGB7ciB4YXJpbmdhbkV4dHJhLWNsaXBib2FyZCwgZWNobyA9IEZBTFNFfQpodG1sdG9vbHM6OnRhZ0xpc3QoCiAgeGFyaW5nYW5FeHRyYTo6dXNlX2NsaXBib2FyZCgKICAgIGJ1dHRvbl90ZXh0ID0gIjxpIGNsYXNzPVwiZmEtc29saWQgZmEtY2xpcGJvYXJkXCIgc3R5bGU9XCJjb2xvcjogIzAwMDA4QlwiPjwvaT4iLAogICAgc3VjY2Vzc190ZXh0ID0gIjxpIGNsYXNzPVwiZmEgZmEtY2hlY2tcIiBzdHlsZT1cImNvbG9yOiAjOTBCRTZEXCI+PC9pPiIsCiAgICBlcnJvcl90ZXh0ID0gIjxpIGNsYXNzPVwiZmEgZmEtdGltZXMtY2lyY2xlXCIgc3R5bGU9XCJjb2xvcjogI0Y5NDE0NFwiPjwvaT4iCiAgKSwKICBybWFya2Rvd246Omh0bWxfZGVwZW5kZW5jeV9mb250X2F3ZXNvbWUoKQopCmBgYAoKCmBgYHtjc3MsIGVjaG8gPSBGQUxTRX0KYm9keSAubWFpbi1jb250YWluZXIgewogIG1heC13aWR0aDogMTAwJSAhaW1wb3J0YW50OwogIHdpZHRoOiAxMDAlICFpbXBvcnRhbnQ7Cn0KYm9keSB7CiAgbWF4LXdpZHRoOiAxMDAlICFpbXBvcnRhbnQ7Cn0KCmJvZHksIHRkIHsKICAgZm9udC1zaXplOiAxNnB4Owp9CmNvZGUucnsKICBmb250LXNpemU6IDE0cHg7Cn0KcHJlIHsKICBmb250LXNpemU6IDE0cHgKfQouY3VzdG9tLWJveCB7CiAgYmFja2dyb3VuZC1jb2xvcjogI2Y1ZjdmYTsgLyogTGlnaHQgZ3JleS1ibHVlIGJhY2tncm91bmQgKi8KICBib3JkZXItY29sb3I6ICNlMWU4ZWQ7IC8qIExpZ2h0IGJvcmRlciBjb2xvciAqLwogIGNvbG9yOiAjMmMzZTUwOyAvKiBEYXJrIHRleHQgY29sb3IgKi8KICBwYWRkaW5nOiAxNXB4OyAvKiBQYWRkaW5nIGluc2lkZSB0aGUgYm94ICovCiAgYm9yZGVyLXJhZGl1czogNXB4OyAvKiBSb3VuZGVkIGNvcm5lcnMgKi8KICBtYXJnaW4tYm90dG9tOiAyMHB4OyAvKiBTcGFjaW5nIGJlbG93IHRoZSBib3ggKi8KfQouY2FwdGlvbiB7CiAgbWFyZ2luOiBhdXRvOwogIHRleHQtYWxpZ246IGNlbnRlcjsKICBtYXJnaW4tYm90dG9tOiAyMHB4OyAvKiBTcGFjaW5nIGJlbG93IHRoZSBib3ggKi8KfQpgYGAKCgpMZXQgdXMgc2V0IHNvbWUgZ2xvYmFsIG9wdGlvbnMgZm9yIGFsbCBjb2RlIGNodW5rcyBpbiB0aGlzIGRvY3VtZW50LgoKCmBgYHtyfQojIFNldCBzZWVkIGZvciByZXByb2R1Y2liaWxpdHkKc2V0LnNlZWQoMTk4MikgCiMgU2V0IGdsb2JhbCBvcHRpb25zIGZvciBhbGwgY29kZSBjaHVua3MKa25pdHI6Om9wdHNfY2h1bmskc2V0KAogICMgRGlzYWJsZSBtZXNzYWdlcyBwcmludGVkIGJ5IFIgY29kZSBjaHVua3MKICBtZXNzYWdlID0gVFJVRSwgICAgCiAgIyBEaXNhYmxlIHdhcm5pbmdzIHByaW50ZWQgYnkgUiBjb2RlIGNodW5rcwogIHdhcm5pbmcgPSBUUlVFLCAgICAKICAjIFNob3cgUiBjb2RlIHdpdGhpbiBjb2RlIGNodW5rcyBpbiBvdXRwdXQKICBlY2hvID0gVFJVRSwgICAgICAgIAogICMgSW5jbHVkZSBib3RoIFIgY29kZSBhbmQgaXRzIHJlc3VsdHMgaW4gb3V0cHV0CiAgaW5jbHVkZSA9IFRSVUUsICAgICAKICAjIEV2YWx1YXRlIFIgY29kZSBjaHVua3MKICBldmFsID0gVFJVRSwgICAgICAgCiAgIyBFbmFibGUgY2FjaGluZyBvZiBSIGNvZGUgY2h1bmtzIGZvciBmYXN0ZXIgcmVuZGVyaW5nCiAgY2FjaGUgPSBGQUxTRSwgICAgICAKICAjIEFsaWduIGZpZ3VyZXMgaW4gdGhlIGNlbnRlciBvZiB0aGUgb3V0cHV0CiAgZmlnLmFsaWduID0gImNlbnRlciIsCiAgIyBFbmFibGUgcmV0aW5hIGRpc3BsYXkgZm9yIGhpZ2gtcmVzb2x1dGlvbiBmaWd1cmVzCiAgcmV0aW5hID0gMiwKICAjIFNob3cgZXJyb3JzIGluIHRoZSBvdXRwdXQgaW5zdGVhZCBvZiBzdG9wcGluZyByZW5kZXJpbmcKICBlcnJvciA9IFRSVUUsCiAgIyBEbyBub3QgY29sbGFwc2UgY29kZSBhbmQgb3V0cHV0IGludG8gYSBzaW5nbGUgYmxvY2sKICBjb2xsYXBzZSA9IEZBTFNFCikKIyBTdGFydCB0aGUgZmlndXJlIGNvdW50ZXIKZmlnX2NvdW50IDwtIDAKIyBEZWZpbmUgdGhlIGNhcHRpb25lciBmdW5jdGlvbgpjYXB0aW9uZXIgPC0gZnVuY3Rpb24oY2FwdGlvbikgewogIGZpZ19jb3VudCA8PC0gZmlnX2NvdW50ICsgMQogIHBhc3RlMCgiRmlndXJlICIsIGZpZ19jb3VudCwgIjogIiwgY2FwdGlvbikKfQojIERlZmluZSB0aGUgZnVuY3Rpb24gdG8gdHJ1bmNhdGUgYSBudW1iZXIgdG8gdHdvIGRlY2ltYWwgcGxhY2VzCnRydW5jYXRlX3RvX3R3byA8LSBmdW5jdGlvbih4KSB7CiAgZmxvb3IoeCAqIDEwMCkgLyAxMDAKfQpgYGAKCmBgYHtyfQpteS5nZXQucm9vdHMgPC0gZnVuY3Rpb24ob3JkZXIsIGJldGEpIHsKICBtdCA8LSBnZXQocGFzdGUwKCJtIiwgb3JkZXIsICJ0YWJsZSIpKQogIHJiIDwtIHJlcCgwLCBvcmRlciArIDEpCiAgcmMgPC0gcmVwKDAsIG9yZGVyKQogIGlmKG9yZGVyID09IDEpIHsKICAgICAgcmMgPSBhcHByb3gobXQkYmV0YSwgbXRbW3Bhc3RlMCgicmMiKV1dLCBiZXRhKSR5CiAgICB9IGVsc2UgewogICAgICByYyA9IHNhcHBseSgxOm9yZGVyLCBmdW5jdGlvbihpKSB7CiAgICAgICAgYXBwcm94KG10JGJldGEsIG10W1twYXN0ZTAoInJjLiIsIGkpXV0sIGJldGEpJHkKICAgICAgfSkKICAgIH0KICAgIHJiID0gc2FwcGx5KDE6KG9yZGVyKzEpLCBmdW5jdGlvbihpKSB7CiAgICAgIGFwcHJveChtdCRiZXRhLCBtdFtbcGFzdGUwKCJyYi4iLCBpKV1dLCB4b3V0ID0gYmV0YSkkeQogICAgfSkKICAgIGZhY3RvciA9IGFwcHJveChtdCRiZXRhLCBtdCRmYWN0b3IsIHhvdXQgPSBiZXRhKSR5CiAgcmV0dXJuKGxpc3QocmIgPSByYiwgcmMgPSByYywgZmFjdG9yID0gZmFjdG9yKSkKfQoKbTF0YWJsZSA8LSByU1BERTo6Om0xdGFibGUKbTJ0YWJsZSA8LSByU1BERTo6Om0ydGFibGUKbTN0YWJsZSA8LSByU1BERTo6Om0zdGFibGUKbTR0YWJsZSA8LSByU1BERTo6Om00dGFibGUKCnBvbHlfZnJvbV9yb290cyA8LSBmdW5jdGlvbihyb290cykgewogIGNvZWYgPC0gMQogIGZvciAociBpbiByb290cykgewogICAgY29lZiA8LSBjb252b2x2ZShjb2VmLCBjKDEsIC1yKSwgdHlwZSA9ICJvcGVuIikKICB9CiAgcmV0dXJuKGNvZWYpICMgcmV0dXJucyBpbiBpbmNyZWFzaW5nIG9yZGVyIGxpa2UgYStieCtjeF4yKy4uLgp9CiAKIyBHcm91cCBjb21wbGV4IHJvb3RzIGludG8gY29uanVnYXRlIHBhaXJzCmdldF9jb25qdWdhdGVfcGFpcnMgPC0gZnVuY3Rpb24ocm9vdHMpIHsKICB1c2VkIDwtIHJlcChGQUxTRSwgbGVuZ3RoKHJvb3RzKSkKICBwYWlycyA8LSBsaXN0KCkKICBmb3IgKGkgaW4gc2VxX2Fsb25nKHJvb3RzKSkgewogICAgaWYgKCF1c2VkW2ldKSB7CiAgICAgIGNvbmpfcm9vdCA8LSBDb25qKHJvb3RzW2ldKQogICAgICAjIEZpbmQgaXRzIGNvbmp1Z2F0ZSAod2l0aGluIHRvbGVyYW5jZSkKICAgICAgaiA8LSB3aGljaCghdXNlZCAmIGFicyhSZShyb290cykgLSBSZShjb25qX3Jvb3QpKSA8IDFlLTggJiBhYnMoSW0ocm9vdHMpIC0gSW0oY29ual9yb290KSkgPCAxZS04KQogICAgICBqIDwtIGpbaiAhPSBpXQogICAgICBpZiAobGVuZ3RoKGopID4gMCkgewogICAgICAgIHVzZWRbYyhpLCBqWzFdKV0gPC0gVFJVRQogICAgICAgIHBhaXJzW1tsZW5ndGgocGFpcnMpICsgMV1dIDwtIGMocm9vdHNbaV0sIHJvb3RzW2pbMV1dKQogICAgICB9CiAgICB9CiAgfQogIHJldHVybihwYWlycykKfQoKY29tcHV0ZV9zdW1fcG9seSA8LSBmdW5jdGlvbihmYWN0b3IsIHByX3Jvb3RzLCBwbF9yb290cywgY3RlKSB7CiAgCiAgcHJfY29lZiA8LSBjKDAsIHBvbHlfZnJvbV9yb290cyhwcl9yb290cykpICMgaW4gZGVjcmVhc2luZyBvcmRlciBsaWtlIHhebitieF4obi0xKStjeF4obi0yKSsuLi4KICBwbF9jb2VmIDwtIHBvbHlfZnJvbV9yb290cyhwbF9yb290cykgIyBpbiBkZWNyZWFzaW5nIG9yZGVyIGxpa2UgeF5uK2J4XihuLTEpK2N4XihuLTIpKy4uLgogIHByX3BsdXNfcGxfY29lZiA8LSBmYWN0b3IgKiBwcl9jb2VmICsgY3RlICogcGxfY29lZgogIHJldHVybihwcl9wbHVzX3BsX2NvZWYpCn0KCmNvbXB1dGVfcmVhbF9yb290c19hbmRfY29tcGxleF9jb2VmIDwtIGZ1bmN0aW9uKHByX3BsdXNfcGxfY29lZikgewoKICBwcl9wbHVzX3BsX3Jvb3RzIDwtIHBvbHlyb290KHJldihwcl9wbHVzX3BsX2NvZWYpKQogIAogIHJlYWxfcm9vdHMgPC0gUmUocHJfcGx1c19wbF9yb290c1thYnMoSW0ocHJfcGx1c19wbF9yb290cykpIDwgMWUtOF0pCiAgY29tcGxleF9yb290cyA8LSBwcl9wbHVzX3BsX3Jvb3RzW2FicyhJbShwcl9wbHVzX3BsX3Jvb3RzKSkgPj0gMWUtOF0KICAKICBjb21wbGV4X3BvbHlfY29lZnMgPC0gbGlzdCgpCiAgaWYgKGxlbmd0aChjb21wbGV4X3Jvb3RzKSA+IDApIHsKICAgIHBhaXJzIDwtIGdldF9jb25qdWdhdGVfcGFpcnMoY29tcGxleF9yb290cykKICAgIGZvciAocGFpciBpbiBwYWlycykgewogICAgICBjb2VmIDwtIHJldihSZShwb2x5X2Zyb21fcm9vdHMoYyhwYWlyWzFdLCBwYWlyWzJdKSkpKSAgIyByZXR1cm5zIGluIHheMiArIGJ4ICsgYyBvcmRlcgogICAgICBjb21wbGV4X3BvbHlfY29lZnNbW2xlbmd0aChjb21wbGV4X3BvbHlfY29lZnMpICsgMV1dIDwtIHJvdW5kKGNvZWYsIDEyKQogICAgfQogIH0gZWxzZSB7CiAgICBjb21wbGV4X3BvbHlfY29lZnMgPC0gbGlzdCgpICAjIE5vIGNvbXBsZXggcm9vdHMKICB9CiAgCiAgcmV0dXJuKGxpc3QobmV3X2ZhY3RvciA9IHByX3BsdXNfcGxfY29lZlsxXSwgCiAgICAgICAgICAgICAgcHJfcGx1c19wbF9yb290cyA9IHByX3BsdXNfcGxfcm9vdHMsIAogICAgICAgICAgICAgIHJlYWxfcm9vdHMgPSByZWFsX3Jvb3RzLCAKICAgICAgICAgICAgICBjb21wbGV4X3BvbHlfY29lZnMgPSBjb21wbGV4X3BvbHlfY29lZnMpKQp9CgpteS5mcmFjdGlvbmFsLm9wZXJhdG9ycyA8LSBmdW5jdGlvbihMLCAjIGthcHBhXjJDICsgRwogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBiZXRhLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBDLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBzY2FsZS5mYWN0b3IsICMga2FwcGFeMgogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtID0gMSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdGltZV9zdGVwKSB7CgogIEMgPC0gTWF0cml4OjpEaWFnb25hbChkaW0oQylbMV0sIHJvd1N1bXMoQykpICMgbHVtcGVkCiAgQ2kgPC0gTWF0cml4OjpEaWFnb25hbChkaW0oQylbMV0sIDEgLyByb3dTdW1zKEMpKSAjIGx1bXBlZCAKICBJIDwtIE1hdHJpeDo6RGlhZ29uYWwoZGltKEMpWzFdKQogIEwgPC0gTCAvIHNjYWxlLmZhY3RvciAjIEMgKyBHL2thcHBhXjIKICBMQ2kgPC0gTCAlKiUgQ2kKICByb290cyA8LSBteS5nZXQucm9vdHMobSwgYmV0YSkKICBQbC5yb290cyA8LSByb290cyRyYgogIFByLnJvb3RzIDwtIHJvb3RzJHJjCiAgZmFjdG9yIDwtIHJvb3RzJGZhY3RvcgogIG5ld19mYWN0b3JfYW5kX3Jvb3RzIDwtIGNvbXB1dGVfcmVhbF9yb290c19hbmRfY29tcGxleF9jb2VmKGNvbXB1dGVfc3VtX3BvbHkoZmFjdG9yLCBQci5yb290cywgUGwucm9vdHMsIHRpbWVfc3RlcCkpCiAgCiAgbmV3X3JlYWxfcm9vdHMgPC0gbmV3X2ZhY3Rvcl9hbmRfcm9vdHMkcmVhbF9yb290cwogIG5ld19mYWN0b3IgPC0gbmV3X2ZhY3Rvcl9hbmRfcm9vdHMkbmV3X2ZhY3RvcgogIGNvbXBsZXhfcG9seV9jb2VmcyA8LSBuZXdfZmFjdG9yX2FuZF9yb290cyRjb21wbGV4X3BvbHlfY29lZnMKICAKICBQcl9wbHVzX1BsLmZhY3RvcnMgPC0gbGlzdCgpCiAgaWYgKGxlbmd0aChuZXdfcmVhbF9yb290cykgPj0gMSkgewogICAgZm9yIChpIGluIDE6bGVuZ3RoKG5ld19yZWFsX3Jvb3RzKSkgewogICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbW2ldXSA8LSBMQ2kgLSBuZXdfcmVhbF9yb290c1tpXSAqIEkKICAgIH0KICB9CiAgbGVuZ3RoX25vdyA8LSBsZW5ndGgoUHJfcGx1c19QbC5mYWN0b3JzKQogIGlmKGxlbmd0aChjb21wbGV4X3BvbHlfY29lZnMpID49IDEpIHsKICAgIGZvciAoaSBpbiAxOmxlbmd0aChjb21wbGV4X3BvbHlfY29lZnMpKSB7CiAgICB9CiAgICBQcl9wbHVzX1BsLmZhY3RvcnNbW2xlbmd0aF9ub3cgKyBpXV0gPC0gTENpICUqJSBMQ2kgKyBjb21wbGV4X3BvbHlfY29lZnNbW2ldXVsyXSAqIExDaSArIGNvbXBsZXhfcG9seV9jb2Vmc1tbaV1dWzNdICogSQogIH0KICAKICBQci5mYWN0b3JzIDwtIGxpc3QoKQogIFByLmZhY3RvcnNbWzFdXSA8LSBJIC0gTENpICogcm9vdHMkcmNbMV0KCiAgaWYgKGxlbmd0aChyb290cyRyYykgPiAxKSB7CiAgICBmb3IgKGkgaW4gMjpsZW5ndGgocm9vdHMkcmMpKSB7CiAgICAgIFByLmZhY3RvcnNbW2ldXSA8LSBJIC0gTENpICogcm9vdHMkcmNbaV0KICAgIH0KICB9CiAgCiAgI1ByLmZhY3RvcnNbW2xlbmd0aChQci5mYWN0b3JzKSArIDFdXSA8LSBJCiAgIAogIG91dHB1dCA8LSBsaXN0KAogICAgQ2kgPSBDaSwKICAgIEMgPSBDLAogICAgTENpID0gTENpLAogICAgTCA9IEwsCiAgICBtID0gbSwKICAgIGJldGEgPSBiZXRhLAogICAgZmFjdG9yID0gZmFjdG9yLAogICAgbmV3X2ZhY3RvciA9IG5ld19mYWN0b3IsCiAgICBQci5mYWN0b3JzID0gUHIuZmFjdG9ycywKICAgIFByX3BsdXNfUGwuZmFjdG9ycyA9IFByX3BsdXNfUGwuZmFjdG9ycwogICkKICByZXR1cm4ob3V0cHV0KQp9CmBgYAoKCgoKYGBge3J9CiMgaW5zdGFsbC5wYWNrYWdlcygiSU5MQSIscmVwb3M9YyhnZXRPcHRpb24oInJlcG9zIiksSU5MQT0iaHR0cHM6Ly9pbmxhLnItaW5sYS1kb3dubG9hZC5vcmcvUi90ZXN0aW5nIiksIGRlcD1UUlVFKQojIGlubGEudXBncmFkZSh0ZXN0aW5nID0gVFJVRSkKIyByZW1vdGVzOjppbnN0YWxsX2dpdGh1YigiaW5sYWJydS1vcmcvaW5sYWJydSIsIHJlZiA9ICJkZXZlbCIpCiMgcmVtb3Rlczo6aW5zdGFsbF9naXRodWIoImRhdmlkYm9saW4vcnNwZGUiLCByZWYgPSAiZGV2ZWwiKQojIHJlbW90ZXM6Omluc3RhbGxfZ2l0aHViKCJkYXZpZGJvbGluL21ldHJpY2dyYXBoIiwgcmVmID0gImRldmVsIikKbGlicmFyeShJTkxBKQpsaWJyYXJ5KGlubGFicnUpCmxpYnJhcnkoclNQREUpCmxpYnJhcnkoTWV0cmljR3JhcGgpCmxpYnJhcnkoZ3JhdGVmdWwpCgpsaWJyYXJ5KHBsb3RseSkKYGBgCgoKV2Ugd2FudCB0byBzb2x2ZSB0aGUgZnJhY3Rpb25hbCBkaWZmdXNpb24gZXF1YXRpb24KXGJlZ2lue2VxdWF0aW9ufQpcbGFiZWx7ZXE6bWFpbmVxfQogICAgXHBhcnRpYWxfdCB1Kyhca2FwcGFeMi1cRGVsdGFfXEdhbW1hKV57XGZyYWN7XGFscGhhfXsyfX0gdT1mIFx0ZXh0IHsgb24gfSBcR2FtbWEgXHRpbWVzKDAsIFQpLCBccXVhZCB1KDApPXVfMCBcdGV4dCB7IG9uIH0gXEdhbW1hLApcZW5ke2VxdWF0aW9ufQp3aGVyZSAkdSQgc2F0aXNmaWVzIHRoZSBLaXJjaGhvZmYgdmVydGV4IGNvbmRpdGlvbnMKXGJlZ2lue2VxdWF0aW9ufQpcbGFiZWx7ZXE6S2NvbmR9CiAgICBcbGVmdFx7XHBoaVxpbiBDKFxHYW1tYSlcO1xCaWd8XDsgXGZvcmFsbCB2XGluIFY6IFxzdW1fe2VcaW5cbWF0aGNhbHtFfV92fVxwYXJ0aWFsX2UgXHBoaSh2KT0wIFxyaWdodFx9ClxlbmR7ZXF1YXRpb259CgpJZiAkZj0wJCwgdGhlbiB0aGUgc29sdXRpb24gaXMgZ2l2ZW4gYnkKXGJlZ2lue2VxdWF0aW9ufQpcbGFiZWx7ZXE6c29sX3JlcHJlbnRhdGlvbn0KICAgICAgICB1KHMsdCkgPSBcZGlzcGxheXN0eWxlXHN1bV97alxpblxtYXRoYmJ7Tn19ZV57LVxsYW1iZGFee1xmcmFje1xhbHBoYX17Mn19X2p0fVxsZWZ0KHVfMCwgZV9qXHJpZ2h0KV97TF8yKFxHYW1tYSl9ZV9qKHMpLgpcZW5ke2VxdWF0aW9ufQoKYGBge3J9CiMgRnVuY3Rpb24gdG8gYnVpbGQgYSB0YWRwb2xlIGdyYXBoIGFuZCBjcmVhdGUgYSBtZXNoCmdldHNfZ3JhcGhfdGFkcG9sZSA8LSBmdW5jdGlvbihoKXsKICBlZGdlMSA8LSByYmluZChjKDAsMCksYygxLDApKQogIHRoZXRhIDwtIHNlcShmcm9tPS1waSx0bz1waSxsZW5ndGgub3V0ID0gMTAwKQogIGVkZ2UyIDwtIGNiaW5kKDErMS9waStjb3ModGhldGEpL3BpLHNpbih0aGV0YSkvcGkpCiAgZWRnZXMgPSBsaXN0KGVkZ2UxLCBlZGdlMikKICBncmFwaCA8LSBtZXRyaWNfZ3JhcGgkbmV3KGVkZ2VzID0gZWRnZXMpCiAgZ3JhcGgkYnVpbGRfbWVzaChoID0gaCkKICByZXR1cm4oZ3JhcGgpCn0KYGBgCgpMZXQgJFxHYW1tYV9UID0gKFxWY2FsLFxFY2FsKSQgY2hhcmFjdGVyaXplIHRoZSB0YWRwb2xlIGdyYXBoIHdpdGggJFxWY2FsID0gXHt2XzEsdl8yXH0kIGFuZCAkXEVjYWwgPSBce2VfMSxlXzJcfSQuIFRoZSBsZWZ0IGVkZ2UgJGVfMSQgaGFzIGxlbmd0aCAxIGFuZCB0aGUgY2lyY3VsYXIgZWRnZSAkZV8yJCBoYXMgbGVuZ3RoIDIuIEFzIGRpc2N1c3NlZCBiZWZvcmUsIGEgcG9pbnQgb24gJGVfMSQgaXMgcGFyYW1ldGVyaXplZCB2aWEgJHM9XGxlZnQoZV8xLCB0XHJpZ2h0KSQgZm9yICR0IFxpblswLDFdJCBhbmQgYSBwb2ludCBvbiAkZV8yJCB2aWEgJHM9XGxlZnQoZV8yLCB0XHJpZ2h0KSQgZm9yICR0XGluWzAsMl0kLiBPbmUgY2FuIHZlcmlmeSB0aGF0ICQtXERlbHRhX1xHYW1tYSQgaGFzIGVpZ2VudmFsdWVzICQwLFxsZWZ0XHsoaSBccGkgLyAyKV4yXHJpZ2h0XH1fe2kgXGluIFxtYXRoYmJ7Tn19JCBhbmQgJFxsZWZ0XHsoaSBccGkgLyAyKV4yXHJpZ2h0XH1fezIgaSBcaW4gXG1hdGhiYntOfX0kIHdpdGggY29ycmVzcG9uZGluZyBlaWdlbmZ1bmN0aW9ucyAkXHBoaV8wJCwgJFxsZWZ0XHtccGhpX2lccmlnaHRcfV97aSBcaW4gXG1hdGhiYntOfX0kLCBhbmQgJFxsZWZ0XHtccHNpX2lccmlnaHRcfV97MiBpIFxpbiBcbWF0aGJie059fSQgZ2l2ZW4gYnkgJFxwaGlfMChzKT0xIC8gXHNxcnR7M30kIGFuZCAKXGJlZ2lue2VxdWF0aW9uKn0KICAgIFxwaGlfaShzKT1DX3tccGhpLCBpfVxiZWdpbntjYXNlc30KICAgICAgICAtMiBcc2luIChcZnJhY3tpXHBpfXsyfSkgXGNvcyAoXGZyYWN7aSBccGkgdH17Mn0pLCAmIHMgXGluIGVfMSwgXFwKXHNpbiAoaSBccGkgdCAvIDIpLCAmIHMgXGluIGVfMiwKICAgIFxlbmR7Y2FzZXN9LApccXVhZCAKICAgIFxwc2lfaShzKT1cZnJhY3tcc3FydHszfX17XHNxcnR7Mn19IFxiZWdpbntjYXNlc30KICAgICgtMSlee2kgLyAyfSBcY29zIChcZnJhY3tpIFxwaSB0fXsyfSksICYgcyBcaW4gZV8xLCBcXApcY29zIChcZnJhY3tpIFxwaSB0fXsyfSksICYgcyBcaW4gZV8yLApcZW5ke2Nhc2VzfSwKXGVuZHtlcXVhdGlvbip9CndoZXJlICRDX3tccGhpLCBpfT0xJCBpZiAkaSQgaXMgZXZlbiBhbmQgJENfe1xwaGksIGl9PTEgLyBcc3FydHszfSQgb3RoZXJ3aXNlLiBNb3Jlb3ZlciwgdGhlc2UgZnVuY3Rpb25zIGZvcm0gYW4gb3J0aG9ub3JtYWwgYmFzaXMgZm9yICRMXzIoXEdhbW1hX1QpJC4KCmBgYHtyfQojIEZ1bmN0aW9uIHRvIGNvbXB1dGUgdGhlIGVpZ2VuZnVuY3Rpb25zIAp0YWRwb2xlLmVpZyA8LSBmdW5jdGlvbihrLGdyYXBoKXsKeDEgPC0gYygwLGdyYXBoJGdldF9lZGdlX2xlbmd0aHMoKVsxXSpncmFwaCRtZXNoJFB0RVtncmFwaCRtZXNoJFB0RVssMV09PTEsMl0pIAp4MiA8LSBjKDAsZ3JhcGgkZ2V0X2VkZ2VfbGVuZ3RocygpWzJdKmdyYXBoJG1lc2gkUHRFW2dyYXBoJG1lc2gkUHRFWywxXT09MiwyXSkgCgppZihrPT0wKXsgCiAgZi5lMSA8LSByZXAoMSxsZW5ndGgoeDEpKSAKICBmLmUyIDwtIHJlcCgxLGxlbmd0aCh4MikpIAogIGYxID0gYyhmLmUxWzFdLGYuZTJbMV0sZi5lMVstMV0sIGYuZTJbLTFdKSAKICBmID0gbGlzdChwaGk9ZjEvc3FydCgzKSkgCiAgCn0gZWxzZSB7CiAgZi5lMSA8LSAtMipzaW4ocGkqayoxLzIpKmNvcyhwaSprKngxLzIpIAogIGYuZTIgPC0gc2luKHBpKmsqeDIvMikgICAgICAgICAgICAgICAgICAKICAKICBmMSA9IGMoZi5lMVsxXSxmLmUyWzFdLGYuZTFbLTFdLCBmLmUyWy0xXSkgCiAgCiAgaWYoKGsgJSUgMik9PTEpeyAKICAgIGYgPSBsaXN0KHBoaT1mMS9zcXJ0KDMpKSAKICB9IGVsc2UgeyAKICAgIGYuZTEgPC0gKC0xKV57ay8yfSpjb3MocGkqayp4MS8yKQogICAgZi5lMiA8LSBjb3MocGkqayp4Mi8yKQogICAgZjIgPSBjKGYuZTFbMV0sZi5lMlsxXSxmLmUxWy0xXSxmLmUyWy0xXSkgCiAgICBmIDwtIGxpc3QocGhpPWYxLHBzaT1mMi9zcXJ0KDMvMikpCiAgfQp9CgpyZXR1cm4oZikKfQpgYGAKCkltcGxlbWVudGF0aW9uIG9mICR1JAoKYGBge3J9CmggPC0gMC4wMQpncmFwaCA8LSBnZXRzX2dyYXBoX3RhZHBvbGUoaCA9IGgpClRfZmluYWwgPC0gMgp0aW1lX3N0ZXAgPC0gMC4wMQp0aW1lX3NlcSA8LSBzZXEoMCwgVF9maW5hbCwgYnkgPSB0aW1lX3N0ZXApCiMgQ29tcHV0ZSB0aGUgRkVNIG1hdHJpY2VzCmdyYXBoJGNvbXB1dGVfZmVtKCkKRyA8LSBncmFwaCRtZXNoJEcKQyA8LSBncmFwaCRtZXNoJEMKSSA8LSBNYXRyaXg6OkRpYWdvbmFsKG5yb3coQykpCnggPC0gZ3JhcGgkbWVzaCRWWywgMV0KeSA8LSBncmFwaCRtZXNoJFZbLCAyXQplZGdlX251bWJlciA8LSBncmFwaCRtZXNoJFZ0RVssIDFdCnBvcyA8LSBzdW0oZWRnZV9udW1iZXIgPT0gMSkrMQpvcmRlcl90b19wbG90IDwtIGZ1bmN0aW9uKHYpcmV0dXJuKGModlsxXSwgdlszOnBvc10sIHZbMl0sIHZbKHBvcysxKTpsZW5ndGgodildLCB2WzJdKSkKd2VpZ2h0cyA8LSBncmFwaCRtZXNoJHdlaWdodHMKIyBJbml0aWFsIGNvbmRpdGlvbgojIFVfMCA8LSAxMCpleHAoLSgoeC0xKV4yICsgKHkpXjIpKQpgYGAKCgpgYGB7cn0Ka2FwcGEgPC0gMQphbHBoYSA8LSAwLjUgIyBmcm9tIDAuNSB0byAyCm0gPSA0CmJldGEgPC0gYWxwaGEvMgpMIDwtIGthcHBhXjIqQyArIEcKYGBgCgoKYGBge3J9Cm15X29wIDwtIG15LmZyYWN0aW9uYWwub3BlcmF0b3JzKEwsIGJldGEsIEMsIHNjYWxlLmZhY3RvciA9IGthcHBhXjIsIG0gPSBtLCB0aW1lX3N0ZXApCgpteS5zb2x2ZXIgPC0gZnVuY3Rpb24ob2JqLCB2KXsKICBQci5mYWN0b3JzIDwtIG9iaiRQci5mYWN0b3JzCiAgUHJfcGx1c19QbC5mYWN0b3JzIDwtIG9iaiRQcl9wbHVzX1BsLmZhY3RvcnMKICBtIDwtIG9iaiRtCiAgQyA8LSBvYmokQwogIENpIDwtIG9iaiRDaQogIGZhY3RvciA8LSBvYmokZmFjdG9yCiAgbmV3X2ZhY3RvciA8LSBvYmokbmV3X2ZhY3RvcgogIGlmIChtPT0xKXsKICAgIHRlbXAgPC0gc29sdmUoUHJfcGx1c19QbC5mYWN0b3JzW1syXV0sIAogICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1sxXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbMV1dLCAKICAgICAgICAgICAgICAgICAgICBDJSoldikpCiAgfSBlbHNlIGlmIChtPT0yKXsKICAgIHRlbXAgPC0gc29sdmUoUHJfcGx1c19QbC5mYWN0b3JzW1szXV0sIAogICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1syXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbMl1dLCAKICAgICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1sxXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgICAgUHJfcGx1c19QbC5mYWN0b3JzW1sxXV0sIAogICAgICAgICAgICAgICAgICAgICAgQyUqJXYpKSkKICB9IGVsc2UgaWYgKG09PTMpewogICAgdGVtcCA8LSBzb2x2ZShQcl9wbHVzX1BsLmZhY3RvcnNbWzRdXSwgCiAgICAgICAgICAgICAgICAgIFByLmZhY3RvcnNbWzNdXSAlKiUgc29sdmUoCiAgICAgICAgICAgICAgICAgICAgUHJfcGx1c19QbC5mYWN0b3JzW1szXV0sIAogICAgICAgICAgICAgICAgICAgIFByLmZhY3RvcnNbWzJdXSAlKiUgc29sdmUoCiAgICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzJdXSwgCiAgICAgICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1sxXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzFdXSwgCiAgICAgICAgICAgICAgICAgICAgICAgIEMlKiV2KSkpKQogIH0gZWxzZSBpZiAobT09NCl7CiAgICB0ZW1wIDwtIHNvbHZlKFByX3BsdXNfUGwuZmFjdG9yc1tbNV1dLCAKICAgICAgICAgICAgICAgICAgUHIuZmFjdG9yc1tbNF1dICUqJSBzb2x2ZSgKICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzRdXSwgCiAgICAgICAgICAgICAgICAgICAgUHIuZmFjdG9yc1tbM11dICUqJSBzb2x2ZSgKICAgICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbM11dLCAKICAgICAgICAgICAgICAgICAgICAgIFByLmZhY3RvcnNbWzJdXSAlKiUgc29sdmUoCiAgICAgICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbMl1dLCAKICAgICAgICAgICAgICAgICAgICAgICAgUHIuZmFjdG9yc1tbMV1dICUqJSBzb2x2ZSgKICAgICAgICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzFdXSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgQyUqJXYpKSkpKQogIH0KICByZXR1cm4oKGZhY3Rvci9uZXdfZmFjdG9yKSAqIENpICUqJSB0ZW1wKQp9CmBgYAoKCmBgYHtyfQpvcCA8LSBmcmFjdGlvbmFsLm9wZXJhdG9ycyhMLCBiZXRhLCBDLCBzY2FsZS5mYWN0b3IgPSBrYXBwYV4yLCBtID0gbSkKUGwgPC0gb3AkUGwKUHIgPC0gb3AkUHIKQ2kgPC0gb3AkQ2kKQyA8LSBvcCRDCmBgYAoKCgpgYGB7cn0KIyBQYXJhbWV0ZXJzIHRvIGNvbnN0cnVjdCBVXzAKTl9maW5pdGUgPC0gNCAjIGNob29zZSBhbiBldmVuIG51bWJlcgphZGp1c3RlZF9OX2Zpbml0ZSA8LSBOX2Zpbml0ZSArIE5fZmluaXRlLzIgKyAxCkVJR0VOVkFMX0FMUEhBIDwtIE5VTEwKRUlHRU5GVU4gPC0gTlVMTCAgICAgICAKZm9yIChqIGluIDA6Tl9maW5pdGUpIHsKICAgIGxhbWJkYV9qX2FscGhhIDwtIChrYXBwYV4yICsgKGoqcGkvMileMileKGFscGhhLzIpCiAgICBlX2ogPC0gdGFkcG9sZS5laWcoaixncmFwaCkkcGhpCiAgICBFSUdFTlZBTF9BTFBIQSA8LSBjKEVJR0VOVkFMX0FMUEhBLCBsYW1iZGFfal9hbHBoYSkgICAgICAgICAKICAgIEVJR0VORlVOIDwtIGNiaW5kKEVJR0VORlVOLCBlX2opICAgICAgICAgICAgCiAgICBpZiAoaj4wICYmIChqICUlIDIgPT0gMCkpIHsKICAgICAgbGFtYmRhX2pfYWxwaGEgPC0gKGthcHBhXjIgKyAoaipwaS8yKV4yKV4oYWxwaGEvMikKICAgICAgZV9qIDwtIHRhZHBvbGUuZWlnKGosZ3JhcGgpJHBzaQogICAgICBFSUdFTlZBTF9BTFBIQSA8LSBjKEVJR0VOVkFMX0FMUEhBLCBsYW1iZGFfal9hbHBoYSkgICAgICAgICAKICAgICAgRUlHRU5GVU4gPC0gY2JpbmQoRUlHRU5GVU4sIGVfaikgICAgICAgICAgIAogICAgfQp9CgojIEJ1aWxkaW5nIHRoZSBpbml0aWFsIGNvbmRpdGlvbiBhcyBcc3VtIGNvZWZmX2ogRUlHRU5GVU5fagpjb2VmZiA8LSAxKigxOmFkanVzdGVkX05fZmluaXRlKV4tMQpjb2VmZlstNV0gPC0gMAojIGxvd2VyX3plcm9lciA8LSAxMCAgIyBjaG9vc2UgYW4gZXZlbiBudW1iZXIKIyBhZGp1c3RlZF9sb3dlcl96ZXJvZXIgPC0gbG93ZXJfemVyb2VyICsgbG93ZXJfemVyb2VyLzIgKyAxCiMgY29lZmZbYWRqdXN0ZWRfbG93ZXJfemVyb2VyOmFkanVzdGVkX05fZmluaXRlXSA8LSAwClVfMCA8LSBFSUdFTkZVTiAlKiUgY29lZmYKCiMgQnVpbGRpbmcgdGhlIHRydWUgc29sdXRpb24gYXMgXHN1bSBjb2VmZl9qIEVJR0VORlVOX2ogZV57LVxsYW1iZGFfal57XGZyYWN7XGFscGhhfXsyfX10fQpVX3RydWUgPC0gbWF0cml4KE5BLCBucm93ID0gbGVuZ3RoKHgpLCBuY29sID0gbGVuZ3RoKHRpbWVfc2VxKSkKZm9yIChrIGluIDE6bGVuZ3RoKHRpbWVfc2VxKSkgewogIGF1eF9rIDwtIHJlcCgwLCBsZW5ndGgoeCkpCiAgZm9yIChqIGluIDE6YWRqdXN0ZWRfTl9maW5pdGUpIHsKICAgIGF1eF9rIDwtIGF1eF9rICsgZXhwKC10aW1lX3NlcVtrXSpFSUdFTlZBTF9BTFBIQVtqXSkqY29lZmZbal0qRUlHRU5GVU5bLCBqXQogIH0KICBVX3RydWVbLCBrXSA8LSBhdXhfawp9CmBgYAoKYGBge3J9ClVfYXBwcm94IDwtIG1hdHJpeChOQSwgbnJvdyA9IG5yb3coQyksIG5jb2wgPSBsZW5ndGgodGltZV9zZXEpKQpVX2FwcHJveFssIDFdIDwtIFVfMAoKIyBUaW1lLXN0ZXBwaW5nIGxvb3AKZm9yIChrIGluIDE6KGxlbmd0aCh0aW1lX3NlcSkgLSAxKSkgewogIFVfYXBwcm94WywgayArIDFdIDwtIGFzLm1hdHJpeChteS5zb2x2ZXIobXlfb3AsIFVfYXBwcm94Wywga10pKQp9CmBgYAoKYGBgCntyfQpQci5hcHBseS5tdWx0IDwtIGZ1bmN0aW9uKHYpe3JldHVybihQci5tdWx0KG9wLCB2KSl9ClBsLmFwcGx5LnNvbHZlIDwtIGZ1bmN0aW9uKHYpe3JldHVybihQbC5zb2x2ZShvcCwgdikpfQpQbGlDIDwtIGFwcGx5KEMsIDIsIFBsLmFwcGx5LnNvbHZlKSAjIFBsQ14tMQojIFByZWNvbXB1dGUgdGhlIExIUzEgbWF0cml4CmF1eCA8LSBhcHBseShQbGlDLCAyLCBQci5hcHBseS5tdWx0KSAjUHJQbF4tMUMKTEhTIDwtIGF1eCArIHRpbWVfc3RlcCAqIE1hdHJpeDo6RGlhZ29uYWwobnJvdyhDKSkgCgojIEluaXRpYWxpemUgVSBtYXRyaXggdG8gc3RvcmUgc29sdXRpb24gYXQgZWFjaCB0aW1lIHN0ZXAKVV9hcHByb3ggPC0gbWF0cml4KE5BLCBucm93ID0gbnJvdyhDKSwgbmNvbCA9IGxlbmd0aCh0aW1lX3NlcSkpClVfYXBwcm94WywgMV0gPC0gVV8wCgojIFRpbWUtc3RlcHBpbmcgbG9vcApmb3IgKGsgaW4gMToobGVuZ3RoKHRpbWVfc2VxKSAtIDEpKSB7CiAgIyBDb21wdXRlIHRoZSByaWdodC1oYW5kIHNpZGUgZm9yIHRoZSBzZWNvbmQgZXF1YXRpb24KICBSSFMgPC0gYXV4ICUqJSBVX2FwcHJveFssIGtdCiAgVV9hcHByb3hbLCBrICsgMV0gPC0gYXMubWF0cml4KHNvbHZlKExIUywgUkhTKSkKfQpgYGAKCldlIGFycml2ZSBhdCB0aGUgc2NoZW1lCgpcYmVnaW57ZXF1YXRpb259CihQX3JeVEMrXHRhdSBQX1xlbGxeVClVXntrKzF9ID0gUF9yXlRDVV5rClxsYWJlbHtlcTpzY2hlbWV9ClxlbmR7ZXF1YXRpb259CndoZXJlClxiZWdpbntlcXVhdGlvbn0KUF9yID0gY19tXHByb2Rfe2k9MX1ebSAoSS1yX3sxaX1DXnstMX1MKVxxdWFkXHRleHR7YW5kfVxxdWFkIFBfXGVsbCA9IGJfe20rMX1DXHByb2Rfe2o9MX1ee20rMX0gKEktcl97Mmp9Q157LTF9TCkKXGVuZHtlcXVhdGlvbn0KT2JzZXJ2ZSB0aGF0IApcYmVnaW57ZXF1YXRpb259ClBfcl5UID0gY19tXHByb2Rfe2k9bX1eMSAoSS1yX3sxaX1MQ157LTF9KVxxdWFkXHRleHR7YW5kfVxxdWFkIFBfXGVsbF5UID0gYl97bSsxfVxwcm9kX3tqPW0rMX1eezF9IChJLXJfezJqfUxDXnstMX0pXGNkb3QgQwpcZW5ke2VxdWF0aW9ufQpSZXBsYWNpbmcgdGhlc2UgdHdvIGluIG91ciBzY2hlbWUgd2UgZ2V0ClxiZWdpbntlcXVhdGlvbn0KXGxlZnQoY19tXHByb2Rfe2k9bX1eMSAoSS1yX3sxaX1MQ157LTF9KStcdGF1IGJfe20rMX1ccHJvZF97aj1tKzF9XnsxfSAoSS1yX3syan1MQ157LTF9KVxyaWdodClDVV57aysxfSA9IGNfbVxwcm9kX3tpPW19XjEgKEktcl97MWl9TENeey0xfSlcY2RvdCBDVV5rClxlbmR7ZXF1YXRpb259CldlIGNhbiBlcXVpdmFsZW50bHkgd3JpdGUgdGhpcyBhcwpcYmVnaW57ZXF1YXRpb259ClxsZWZ0KFxkZnJhY3tjX219e2Jfe20rMX19XHByb2Rfe2k9MX1ebSAoSS1yX3sxaX1MQ157LTF9KStcdGF1IFxwcm9kX3tqPTF9XnttKzF9IChJLXJfezJqfUxDXnstMX0pXHJpZ2h0KUNVXntrKzF9ID0gXGRmcmFje2NfbX17Yl97bSsxfX1ccHJvZF97aT0xfV5tIChJLXJfezFpfUxDXnstMX0pXGNkb3QgQ1VeawpcZW5ke2VxdWF0aW9ufQoKV2UgbmVlZCB0byBjb21wdXRlIHRoZSByb290cyBvZiB0aGUgcG9seW5vbWlhbApcYmVnaW57ZXF1YXRpb259ClIoeCkgPSBcZGZyYWN7Y19tfXtiX3ttKzF9fSBccHJvZF97aT0xfV5tICh4LXJfezFpfSkrXHRhdSBccHJvZF97aj0xfV57bSsxfSAoeC1yX3syan0pIApcZW5ke2VxdWF0aW9ufQpTYXkgJFIkIGhhcyBsZWFkaW5nIGNvZWZmaWNpZW50ICRcdGF1JCBhbmQgJG0rMSQgcm9vdHMgJFJfayQgZm9yICRrPTEsXGxkb3RzLG0rMSQuIFRoYXQgaXMsClxiZWdpbntlcXVhdGlvbn0KUih4KSA9IFx0YXVccHJvZF97az0xfV57bSsxfSAoeC1SX2spClxlbmR7ZXF1YXRpb259CldlIGNhbiB0aGVuIHdyaXRlIG91ciBzY2hlbWUgYXMgClxiZWdpbntlcXVhdGlvbn0KXGxlZnQoXHRhdVxwcm9kX3trPTF9XnttKzF9IChJLVJfa0xDXnstMX0pXHJpZ2h0KUNVXntrKzF9ID0gXGRmcmFje2NfbX17Yl97bSsxfX1ccHJvZF97aT0xfV5tIChJLXJfezFpfUxDXnstMX0pXGNkb3QgQ1VeawpcZW5ke2VxdWF0aW9ufQpUaGF0IGlzLApcYmVnaW57ZXF1YXRpb259ClVee2srMX0gPSBcZGZyYWN7Y19tfXtiX3ttKzF9fVxkZnJhY3sxfXtcdGF1fUNeey0xfVxwcm9kX3trPTF9XnttKzF9IChJLVJfa0xDXnstMX0pXnstMX1ccHJvZF97aT0xfV5tIChJLXJfezFpfUxDXnstMX0pXGNkb3QgQ1VeawpcZW5ke2VxdWF0aW9ufQoKCldoYXQgaWYgClxiZWdpbntlcXVhdGlvbn0KUih4KSA9IFx0YXUoeF4yK2F4K2IpXHByb2Rfe2s9MX1ee20tMX0gKHgtUl9rKQpcZW5ke2VxdWF0aW9ufQp0aGVuCgpgYGAKe3J9CkxIUyA8LSB0KFByKSUqJSBDICsgdGltZV9zdGVwICogdChQbCkKIyBJbml0aWFsaXplIFUgbWF0cml4IHRvIHN0b3JlIHNvbHV0aW9uIGF0IGVhY2ggdGltZSBzdGVwClVfYXBwcm94IDwtIG1hdHJpeChOQSwgbnJvdyA9IG5yb3coQyksIG5jb2wgPSBsZW5ndGgodGltZV9zZXEpKQpVX2FwcHJveFssIDFdIDwtIFVfMAoKIyBUaW1lLXN0ZXBwaW5nIGxvb3AKZm9yIChrIGluIDE6KGxlbmd0aCh0aW1lX3NlcSkgLSAxKSkgewogICMgQ29tcHV0ZSB0aGUgcmlnaHQtaGFuZCBzaWRlIGZvciB0aGUgc2Vjb25kIGVxdWF0aW9uCiAgUkhTIDwtIHQoUHIpJSolIEMgJSolIFVfYXBwcm94Wywga10KICBVX2FwcHJveFssIGsgKyAxXSA8LSBhcy5tYXRyaXgoc29sdmUoTEhTLCBSSFMpKQp9CmBgYAoKCmBgYHtyfQp4IDwtIG9yZGVyX3RvX3Bsb3QoeCkKeSA8LSBvcmRlcl90b19wbG90KHkpCm1heF9lcnJvcl9hdF9lYWNoX3RpbWUgPC0gYXBwbHkoKFVfdHJ1ZSAtIFVfYXBwcm94KV4yLCAyLCBtZWFuKQoKVV90cnVlIDwtIGFwcGx5KFVfdHJ1ZSwgMiwgb3JkZXJfdG9fcGxvdCkKVV9hcHByb3ggPC0gYXBwbHkoVV9hcHByb3gsIDIsIG9yZGVyX3RvX3Bsb3QpCgoKCiMgQ3JlYXRlIGludGVyYWN0aXZlIHBsb3QKZmlnIDwtIHBsb3RfbHkoKQoKCiMgQWRkIHNlY29uZCBsaW5lIChtYXhfZXJyb3JfYXRfZWFjaF90aW1lKQpmaWcgPC0gZmlnICU+JSBhZGRfdHJhY2UoCiAgeCA9IH50aW1lX3NlcSwgeSA9IH5tYXhfZXJyb3JfYXRfZWFjaF90aW1lLCB0eXBlID0gJ3NjYXR0ZXInLCBtb2RlID0gJ2xpbmVzK21hcmtlcnMnLAogIGxpbmUgPSBsaXN0KGNvbG9yID0gJ3JlZCcsIHdpZHRoID0gMiksCiAgbWFya2VyID0gbGlzdChzaXplID0gNCksCiAgbmFtZSA9ICJNYXggRXJyb3IgVHJ1ZSBhbmQgQXBwcm94IDIiCikKCiMgTGF5b3V0CmZpZyA8LSBmaWcgJT4lIGxheW91dCgKICB0aXRsZSA9ICJNYXggRXJyb3IgYXQgRWFjaCBUaW1lIFN0ZXAiLAogIHhheGlzID0gbGlzdCh0aXRsZSA9ICJUaW1lIiksCiAgeWF4aXMgPSBsaXN0KHRpdGxlID0gIk1heCBFcnJvciIpLAogIGxlZ2VuZCA9IGxpc3QoeCA9IDAuMSwgeSA9IDAuOSkKKQoKCnBsb3RfZGF0YSA8LSBkYXRhLmZyYW1lKAogIHggPSByZXAoeCwgdGltZXMgPSBuY29sKFVfdHJ1ZSkpLAogIHkgPSByZXAoeSwgdGltZXMgPSBuY29sKFVfdHJ1ZSkpLAogIHpfdHJ1ZSA9IGFzLnZlY3RvcihVX3RydWUpLAogIHpfYXBwcm94ID0gYXMudmVjdG9yKFVfYXBwcm94KSwKICBmcmFtZSA9IHJlcCh0aW1lX3NlcSwgZWFjaCA9IGxlbmd0aCh4KSkKKQoKIyBDb21wdXRlIGF4aXMgbGltaXRzCnhfcmFuZ2UgPC0gcmFuZ2UoeCkKeV9yYW5nZSA8LSByYW5nZSh5KQp6X3JhbmdlIDwtIHJhbmdlKGMoVV90cnVlLCBVX2FwcHJveCkpCgojIEluaXRpYWwgcGxvdCBzZXR1cCAoZmlyc3QgZnJhbWUgb25seSkKcCA8LSBwbG90X2x5KHBsb3RfZGF0YSwgZnJhbWUgPSB+ZnJhbWUpICU+JQogIGFkZF90cmFjZSgKICAgIHggPSB+eCwgeSA9IH55LCB6ID0gfnpfdHJ1ZSwKICAgIHR5cGUgPSAic2NhdHRlcjNkIiwgbW9kZSA9ICJsaW5lcyIsCiAgICBuYW1lID0gIlRydWUiLAogICAgbGluZSA9IGxpc3QoY29sb3IgPSAiYmx1ZSIsIHdpZHRoID0gMikKICApICU+JQogIGFkZF90cmFjZSgKICAgIHggPSB+eCwgeSA9IH55LCB6ID0gfnpfYXBwcm94LAogICAgdHlwZSA9ICJzY2F0dGVyM2QiLCBtb2RlID0gImxpbmVzIiwKICAgIG5hbWUgPSAiQXBwcm94IiwKICAgIGxpbmUgPSBsaXN0KGNvbG9yID0gInJlZCIsIHdpZHRoID0gMikKICApICU+JQogIGxheW91dCgKICAgIHNjZW5lID0gbGlzdCgKICAgICAgeGF4aXMgPSBsaXN0KHRpdGxlID0gIngiLCByYW5nZSA9IHhfcmFuZ2UpLAogICAgICB5YXhpcyA9IGxpc3QodGl0bGUgPSAieSIsIHJhbmdlID0geV9yYW5nZSksCiAgICAgIHpheGlzID0gbGlzdCh0aXRsZSA9ICJWYWx1ZSIsIHJhbmdlID0gel9yYW5nZSksCiAgICAgIGFzcGVjdHJhdGlvID0gbGlzdCh4ID0gMi40LCB5ID0gMS4yLCB6ID0gMS4yKSwKICAgICAgICAgICBjYW1lcmEgPSBsaXN0KAogICAgICBleWUgPSBsaXN0KHggPSAxLjUsIHkgPSAxLjUsIHogPSAxKSwgICMgQWRqdXN0IHRoZSB2aWV3cG9pbnQKICAgICAgY2VudGVyID0gbGlzdCh4ID0gMCwgeSA9IDAsIHogPSAwKSkKICAgICksCiAgICB1cGRhdGVtZW51cyA9IGxpc3QoCiAgICAgIGxpc3QoCiAgICAgICAgdHlwZSA9ICJidXR0b25zIiwgc2hvd2FjdGl2ZSA9IEZBTFNFLAogICAgICAgIGJ1dHRvbnMgPSBsaXN0KAogICAgICAgICAgbGlzdChsYWJlbCA9ICJQbGF5IiwgbWV0aG9kID0gImFuaW1hdGUiLAogICAgICAgICAgICAgICBhcmdzID0gbGlzdChOVUxMLCBsaXN0KGZyYW1lID0gbGlzdChkdXJhdGlvbiA9IDEwMCwgcmVkcmF3ID0gVFJVRSksIGZyb21jdXJyZW50ID0gVFJVRSkpKSwKICAgICAgICAgIGxpc3QobGFiZWwgPSAiUGF1c2UiLCBtZXRob2QgPSAiYW5pbWF0ZSIsCiAgICAgICAgICAgICAgIGFyZ3MgPSBsaXN0KE5VTEwsIGxpc3QobW9kZSA9ICJpbW1lZGlhdGUiLCBmcmFtZSA9IGxpc3QoZHVyYXRpb24gPSAwKSwgcmVkcmF3ID0gRkFMU0UpKSkKICAgICAgICApCiAgICAgICkKICAgICksCiAgICB0aXRsZSA9ICJUaW1lOiAwIgogICkKCiMgQ29udmVydCB0byBwbG90bHkgb2JqZWN0IHdpdGggZnJhbWUgaW5mbwpwYiA8LSBwbG90bHlfYnVpbGQocCkKCiMgSW5qZWN0IGN1c3RvbSB0aXRsZXMgaW50byBlYWNoIGZyYW1lCmZvciAoaSBpbiBzZXFfYWxvbmcocGIkeCRmcmFtZXMpKSB7CiAgdCA8LSB0aW1lX3NlcVtpXQogIGVyciA8LSBzaWduaWYobWF4X2Vycm9yX2F0X2VhY2hfdGltZVtpXSwgNCkKICBwYiR4JGZyYW1lc1tbaV1dJGxheW91dCA8LSBsaXN0KHRpdGxlID0gcGFzdGUwKCJUaW1lOiAiLCB0LCAiIHwgTWF4IEVycm9yOiAiLCBlcnIpKQp9CmBgYAoKCmBgYHtyfQpmaWcgICMgRGlzcGxheSB0aGUgcGxvdApgYGAKCgpgYGB7ciwgZmlnLmhlaWdodCA9IDYsIG91dC53aWR0aCA9ICIxMDAlIiwgZmlnLmNhcCA9IGNhcHRpb25lcigiQ2FwdGlvbiIpfQpwYgpgYGAKCgo=